Network Inference Tutorial¶
In this tutorial, we will go through how to use the modelling and inference framework set up in this package.
The first thing we will need to do is importnt the relevant packages:
import netwin as nw
from netwin.models import NetworkFKPP
from netwin import Model, VBModel, infer
import numpy as np
import matplotlib.pyplot as plt
Forward Model¶
The first step is to define a model using the Model abstract base class. This allows the model to interface appropriately with the inference classes. We can do this in the following way:
class MyNetworkFKPPModel(Model):
def f(self, p, t, theta):
k, a = theta
du = k * (-self.L() @ p) + (a * p) * (1 - p)
return du
def solve(self, p, theta):
return odeint(self.f, p, self.t, args=(theta,))
def forward(self, u0):
p = np.exp(u0[:-2])
theta = np.exp(u0[-2:])
u = self.solve(p, theta)
return u
To use this, we need to define an adjacency matrix for the model, since this is needed for the abstract base class to initialise. We can construct a random matrix with numpy.
A = np.reshape(np.random.normal(5, 2, 25), (5, 5))
np.fill_diagonal(A, 0)
A = A / np.max(A)
and then construct the model, here using a predefined version of the network FKPP, like so:
m = NetworkFKPP(A)
We can use this model to run simulations. To do this, we need to define the appropriate parameters.
# random values for initial node parameters
p = np.random.uniform(0,0.5,5)
# set the model parameters along with node parameters
u0 = np.append(p, [2.0, 3.0])
# define the time over which we wish to simulate
m.t = np.linspace(0,1,50)
We can now simulate and plot the model.
sim = m.forward(np.log(u0))
plt.plot(m.t, sim)
plt.show()
Inference¶
Using this simulation as synthetic data, we now wish to infer the values of our model parameters and initial node values.
We are going to do this using variational Bayes (VB). To do this, we set up a VB Model.
# set guess for initial node concentration
p0 = np.ones([5])
# append initial guess for model parameters
u_guess = np.append(p0, [1.0,1.0])
# Instantiate probablistic model of type VBModel
ProbModel = VBModel(model=m, data=sim, init_means=u_guess)
We can then call the generic infer function to optimise our probablistic model.
pm = infer(ProbModel, n=20)
Iteration 0
Iteration 1
Iteration 2
Iteration 3
Iteration 4
Iteration 5
Iteration 6
Iteration 7
Iteration 8
Iteration 9
Iteration 10
Iteration 11
Iteration 12
Iteration 13
Iteration 14
Iteration 15
Iteration 16
Iteration 17
Iteration 18
Iteration 19
Finished!
Plotting the results, we can see that we obtain the correct values for our model parameters:
from netwin import plot_2dmvn
plot_2dmvn(means=np.exp(pm.m()[-2:]), cov=pm.Cov()[-2:,-2:], n=100000)
from netwin import barplot_concentrations
barplot_concentrations(true=p, inferred=np.exp(pm.m()[:-2]))